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Abstract 



Using a realistic model, the mode Griineisen parameters 7 and the temper- 
ature dependent coefficient of linear thermal expansion a(T), are calculated 
for amorphous silicon. The resulting 7 values differ from the crystalline case 
in having all diversity suppressed, except for a minority of high-frequency lo- 
calized and low-frequency resonant modes. The latter have very large, mostly 
negative 7 (up to —31), caused by volume driven internal strain. As a result, 
the values for at(T) are lower than those of crystalline silicon and are sample 
dependent. 
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Unlike thermal conductivity and specific heat, which have "universal" features, thermal 
expansion is a special property of each glass. The coefficient of linear expansion a can be 
either positive or negative, with magnitude sensitive to sample preparation methods [|TJ. 
This is true even at very low temperatures (T <1 K), where a is believed to be associated 
with tunneling modes and large dispersion of 7 values is found 0. 

Here we present an analysis of ot{T) for an atomistic model of amorphous silicon, an 
important electronic material. The value of a(T) has been measured only at T=383 K 

making theoretical modelling particularly useful. Our calculation shows that the value 
of a(T) is lower than that of crystalline silicon. Deviations from the crystalline values 
are more dramatic at lower temperatures, and are caused by volume driven internal strain, 
which makes the thermal expansion sample dependent. The predicted Griineisen parameters 
display a surprising simplicity, which we interpret as evidence for special, possibly generic, 
properties of vibrations in glasses. Specifically, we find that the majority of modes which are 
neither localized nor ballistically propagating (we call them "diff usons" ) , have a property 
of global indistinguishability, whereas high-frequency localized and low-frequency resonant 
modes are distinguished by the special structural imperfections at which they have largest 
amplitudes. In particular, the resonant modes have a considerable dispersion of 7 values 
(from -31 to 4), consistent with our finding that internal strain is largest at the centers of 
these modes. Below IK, very large average magnitudes of 7 (I7I > 10) were measured in 
many glasses ("normal" values are about one). Although the size of our model does not 
allow us to go as low as IK, we show that the resonant modes provide a likely origin of this 
anomaly. Finally, we establish a close link between the mode Griineisen parameters and the 
mode bond-stretching character. 

Insofar as a complete set of normal modes can be defined (i.e., the modes decay on 
timescales much larger than their periods), one expects vibrational entropy S to be well 
approximated by S = k B E a [( n a+1) log(n a +l)-n a logn a ], where n a = (exp(hu a / ksT)-!)^ 1 
is the average equilibrium occupation number of the mode of frequency u a , a being a counting 
label going from 1 to 3N. N is the number of atoms and ks is the Boltzmann constant. The 
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frequency u a (V,T) may depend on both volume V and temperature. Using the standard 
thermodynamic relation 3a = K^^dS '/ '8V)t, with isothermal compressibility k t , we get (see, 
e.g., Ref. §) 

a ( T ) = 7^Y. C *la- (1) 

Here c a = k B (hu a /k B T) 2 n a (n a + l) is the specific heat of a harmonic oscillator, and j a , given 
by 7a = — <91ogu; a /<91og V, are known as the "mode Griineisen parameters"; they measure 
how sensitive the vibrational eigenfrequencies are to the change of volume. Their knowledge 
is essential not only for a(T), but also for the interpretation of the internal friction and 
sound attenuation experiments |||7j. 

Perturbation theory (see, e.g., Ref. ||) gives the following formula for j a , 



1 p a* p a 

(2) 
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Here $^ (^87) are ^ e coefficients of the quadratic (cubic) terms of the Taylor expansion 
of potential energy in terms of the displacements Ui a of atoms i from equilibrium in the 
direction a = (x,y,z), and r ia is the position vector of the zth atom. The vibrational 
eigenvectors e" Q are normalized (J2ia \ e i a \ 2= ^-)i an d Mi is the mass of atom i. The sums in 
Eq. |2| and § are over all atoms. Because of periodicity, in a crystal the label a can be written 
as (Q,A), denoting wave vector and polarization. The eigenvectors e ia (Q, A) are proportional 
to e ?c *' ri / \/N \ and the resulting crystalline phase coherence allows a simplification of the sum 
in Eq. |2]to the neighborhood of a single small unit cell. By contrast, normal modes in a glass 
have no a priori quantum numbers, and quantitative insight is best achieved by numerical 
diagonalization for finite models. 

If a solid is subject to an infinitesimal homogeneous isotropic strain e, its atomic coordi- 
nates change to r' ia = r ia {l + e), and volume to V = V(l + 3e). Unless atoms in the solid 
are fixed by symmetry, r' ia are not the equilibrium coordinates of the solid with volume V. 



This happens in glasses. After the strain is applied, the atoms relax to new equilibrium 
positions r' ia + ef ia . This is the origin of the internal strain (Eq. |3[) occuring during the 
thermal expansion of glasses (and crystals with atoms not in centers of symmetry In 
silicon crystal fi a = 0. 

The interactions between silicon atoms are represented by the Stillinger- Weber (SW) 
potential ||, which performs well when applied to elastic and vibrational properties of 
different silicon phases [|l0| -[l^| . For the atomic coordinates of amorphous silicon we take the 
model introduced and described in detail in Ref. [11|. The algorithm by Wooten, Winer, and 



Weaire (WWW) [T3|] creates a random network structure of silicon atoms, which are further 
relaxed to a local minimum of the SW potential. Here we use 1000 atoms arranged in a 
cube of side 27.543 A, with periodic boundary conditions. The energy/atom of this model 
is —4.102 eV, which is ~ 5% higher than the energy/atom of crystalline silicon with SW 
potential. Also, the model agrees very well with the neutron structure factor S(Q) measured 
by Kugler et al. |TJ] . 

To calculate a(T) we need the values for k^. The corresponding perturbation formulas 
for kt, which take into account internal strain, can be found in Ref. |T3|. We obtain the value 



1.13 x 10~ 12 cm 2 /dyn for the present model of amorphous silicon. For crystalline silicon, 
the SW potential gives the value of 0.986 x 10~ 12 ||16| , agreeing well with the experimental 
value of 1.02 x 10~ 12 [161. Also, for silicon crystal, we must sum over the Q points, which 



we do using the tetrahedron method with 1772 tetrahedra in the irreducible wedge of the 
Brillouin zone. The results are shown in Fig. [1], together with the measured data. 

Compared with experiment fli~7f , a(T) of silicon crystal is reproduced very well at T >200 



K. At lower T the SW potential does not reproduce the observed negative values of a. 
Negative a has been succesfully explained [[L8|] by negative values of 7 for the low energy 
transverse acoustic (TA) branch. As shown in Fig. 0(A), our 7 values are too weakly negative 
in this regime. In Fig. [I] it is predicted that a(T) for amorphous silicon is somewhat lower 
than a(T) for the crystal. For comparison, Fig. [IJ also shows a(T) without considering 
internal strain (fi Q in Eq. |]is set to zero), and the calculation based on a 216 atom model 



of amorphous silicon. Internal strain reduces the values of a by almost 30% at high T. At 
lower T the values become negative. The model dependence is also clear. The 216 atom 
model is based on the same WWW algorithm and SW potentials, but is more topologically 
constrained (e.g., no four-fold rings are allowed), and has lower energy/atom (by w 0.5%) 
and higher density (by ~ 3%) than the 1000 atom model. Its ct(T) is higher and closer to 
the crystalline case and the measured value. This is not surprising. Measurements on silica 
|l|Jl9| showed that annealing history (or density) markedly changed its thermal expansion 
at all temperatures. For example, pure SiC>2 aged at 1400°C and quenched, has ot{T) at 
high T lower by up to 50% than that slowly aged at 1000°C Our calculation predicts that 
future experiments should see similar behavior in amorphous silicon as well. 

To understand the behavior of ot(T) one has to look at the frequency dependence of 
Griineisen parameters. The values of 7 for silicon crystal, for more than 1000 randomly 
chosen Q points from the irreducible wedge of the Brillouin zone, are in Fig. 0(A). Since 
there is a degenerate surface in Q space for a given u, j(Q, A) has a distribution of values at 
each u. These values are further split according to the branch A of corresponding phonons. 
Particularly striking polarization effects appear at low u, where TA phonons form a distinct 
broad band of j a . 

It is instructive to compare Fig. |2|(A) with Fig. ||](B), where we plot the mode bond- 



stretching parameter S a , 
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introduced in Ref. |Mj (we modified slightly the original formula which did not yield the 
same S a for different but degenerate a). The summation is over all atom pairs < i,j >, 
where j is within 2.8 A(the distance where the first nearest neighbor peak in the pair 



correlation function of amorphous silicon ends [fR]]) from i, and is the unit vector in the 
direction of the bond < i,j >. When S a is close to one, mode a is mostly bond-stretching, 
while values closer to zero indicate bond-bending modes. The similarity between the two 
figures is striking. Except at very low u, one is tempted to write "7 a ~ const x S a ." The 



reason for such close relation is that anharmonicity, reflected in j a , is much greater in the 
bond-stretching part of the interatomic interaction than in the bond-bending part. 

We now discuss the behavior of 7 a in amorphous silicon, whose vibrational states were 
discussed earlier JTTIJT^]. At low frequencies (^5 15 meV) the modes are propagons, acoustic- 



phonon-like vibrations propagating ballistically for distances greater than their wavelength. 
Some of the propagons are resonantly trapped at certain places in the sample, with reduced 
amplitude elsewhere (but not exponential decay as in a localized state); but at high u 
there is no continuum of extended states to decay into, and vibrational amplitudes decay 
exponentially with distance. Properties of these resonant modes were studied in different 



glassy systems ||21|| , and for our model will be reported elsewhere |[22|| . Modes between 15 
meV and 71 meV are diffusons, extended states carrying energy diffusively. Since they 
form the majority of the spectrum, diffusons dominate the thermal properties of amorphous 
silicon at temperatures from several kelvins up to the melting point. Above u c ~ 71 meV, 
the mobility edge, the modes are locons, localized states. 

The calculated 7 a are shown in Fig. |3|(A). Three regions can be clearly distinguished, (i) 
At low frequencies the values are very scattered, unusually large, and mostly negative, (ii) 
In the region of diffusons, after a monotonous increase, 7 a becomes almost constant at ^30 
meV, and its values collapse into a very narrow region, (iii) At locon frequencies 7 a spreads 
out again, but much less than at low u. 

Why are the low frequency values of 7 so large? At most one would expect the data 
to be scattered between and 1 due to a partial presence of Q and A for propagons. We 
find that the cause of this anomalous behavior is the resonant modes. As already shown 



2T| , these modes have very large amplitudes at groups of under-coordinated atoms. Similar 
groups of over-coordinated atoms are responsible for locons. Since under-coordination im- 
plies "softness" , it is natural that the internal strain r ia is also largest at these sites. This 
is why the magnitudes of 7 are large for the resonant modes. The dispersion of values, also 
seen for locons, is explained by different nature of topological defects where the modes have 
largest amplitudes. 
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The bands seen in 7(Q, A), Fig. ||(A), of crystalline silicon, are suppressed in the "dif- 
fuson" portion of the spectrum of amorphous silicon. This has the following interpretation. 
In a crystal, knowledge of the pattern of a normal mode in a few adjacent unit cells allows 
one to predict the pattern in distant regions. In a glass, such a prediction is not possible 
without a very large computer calculation. More surprising, different normal modes in the 
"diffuson" part of the spectrum are not globally distinguishable. If one knows the pattern 
of displacement of two modes (al and a2) of similar frequencies in one spatial region, one 
could not recognize which pattern was al and which al in a distant region. There seem not 
to be properties other than u a which can classify these modes. 

Finally, Fig. |3|(B) shows S a for amorphous silicon. Note the surprising uniformity of 
bond-stretching character (see also Ref. p0[), even in the region of propagons where memory 
of polarization character might have been expected. As in the crystalline case, S a closely 
follows 7 a , except at low frequencies where the values of 7 are much more influenced by 
dynamics than by vibrational pattern. 

We thank J. L. Feldman and S. Bickham for stimulating discussions. This work was 
supported by NSF Grant No. DMR 9417755. 
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FIG. 1. Calculated and measured linear thermal expansion a(T) for crystalline and amorphous 
silicon. IS: internal strain. 



FIG. 2. Calculated Griineisen parameters j a (A) and bond stretching parameters S a (B) for 
crystalline silicon as a function of frequency. In (A) the polarization labels are: TA - transverse 
acoustic, LA - longitudinal acoustic, LO - longitudinal optical, TO - transverse optical. The solid 
line is the vibrational densities of states (DOS) in arbitrary units. 



FIG. 3. Calculated Griineisen parameters j a (A) and bond-stretching character S a (B) of 
amorphous silicon as a function of frequency. The vertical line at lo c ps 71 meV is the mobility 
edge. The scale in (A) is split by the horizontal line at —1 to emphasize the large negative values 
of 7 at low to. The solid line is DOS in arbitrary units. 
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